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Abstract 

The aim of this paper is to provide a better understanding of the derivation 
of two moments approximate models of shallow water type for thin power law 
films down an incline. We limit ourselves to the case of laminar flow for which 
the boundary layer issued from the interaction of the flow with the bottom 
surface has an influence all over the transverse direction to the flow. In this 
case the concept itself of thin film and its relation with long wave asymptotic 
leads naturally to flow conditions around a uniform free surface Poiseuille flow. 

This Poiseuille flow exhibits a divergence of the apparent viscosity at the 
free surface which, in turn, introduces a singularity in the formulation of the 
Orr-Sommerfeld equations and in the derivation of shallow water models. We 
remove this singularity by inverting the relation between the fluid strain and 
the deformation tensor and by keeping the fluid strain as a variable. No regu- 
larization procedure is needed here, in contrast to ptQ]. 

In order to derive the shallow water equations, we compute a second order 
asymptotic expansion of the fluid velocity field and fluid strain in the shallow 
water regime through a nonlinear iterative scheme. We show that a linearized 
version of this scheme naturally provides an expansion of the eigenfunctions 
of the Orr-Sommerfeld equations in the small wavenumber /small frequency 
regime. The shallow water models account for the streamwise diffusion of mo- 
mentum which is important to describe accurately the full dynamic of the thin 
film flows when instabilities like roll-waves arise. 
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1 Introduction 

Mathematical models and numerical simulations for the flow of a relatively thin layer 
of Non- Newtonian fluids over an inclined surface and under the effect of gravity have 
important applications in natural processes such as mud flows in mountainous areas 
or submarine continental slopes, debris avalanches: it is of interest to be able to 
predict the characteristics (velocity, pressure, runout extent) of such mud flows for 
geophysical and engineering purposes. 

Herein, we consider the simple situation of a thin liquid layer over an inclined 
plane. This is not a restriction if, instead, one considers an arbitrary topography 
provided that the curvature is asymptotically small and slope never vanishes: in this 
case, we just substitute the constant slope by local slope of the bottom (see the 
discussion in I BCN V] ) . Moreover, since instabilities in such flows are primarily two 
dimensional, we will assume that the flow is two dimensional. Herein, we consider 
an incompressible and Non-Newtonian fluid, modeled by the Ostwald-de Wael (also 
known as the Norton-Hoff law in metallurgy) constitutive law 



li p \\D(u) 



in— 1 



D(u). 



where r is the strain rate of the fluid, u its velocity and D(u) = Vu + Vu T the 
deformation rate. The norm on tensors ||^4|| 2 = Tr(AA ) is the classical Frobenius 
norm. The case n < 1 corresponds to shear thinning fluids (like ice) whereas n > 1 
corresponds to shear thickening fluids (like some polymeric fluids). 



Given the fact that the two dimensional Cauchy Momentum equations are difficult 
to treat both analytically and numerically, in particular if the boundary is free, it is 
important to obtain reduced models that are able to capture the relevant features, 
but are mathematically more manageable. Modeling the motion of a thin films down 
an inclined plane leads to a hierarchy of models. The first approximation level are 
the lubrication models: in the shallow water regime and if inertia is small (small 
Reynold number), the fluid velocity and pressure are determined by the local fluid 
height and its derivatives. In this case one obtains a model for the local fluid height 
in the form d t h + d x G(h, d x h, , d™h) = where G involves various algebraic powers 
and differentiation orders of h. This type of equations is obtained by means of asymp- 
totic expansions in power of e, usually called film parameter or aspect ratio between 
the width of the layer and the typical length of the phenomena in the streamwise 



direction. The simplest models obtained for Newtonian fluids are the Benney's 
and generalized Kuramoto-Sivashinsky equations. The simplification brought by this 
reduction has permitted ft^rst study of the nonlinear development of waves using 
dynamical system theory |PMP]. In the case of Non- Newtonian fluids described by a 
power-law fluids (0 < n < 1 for mud flows or ice sheets dynamic), Ng and Mei [14] 
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obtained a Benney type equation and found a stability criterion for Poiseuille flows 
under large wavelength perturbations. However, Benney's equation exhibits finite 
time singularities and is only valid for small amplitude waves. 

In order to obtain models valid for thicker flows and larger Reynolds number, one 
considers shallow water type flows: a specific shape for the velocity profile is assumed 
together with a hydrostatic field assumption. Averaging the streamwise momentum 
equation and the divergence free condition, one obtains an evolution system for the 
local fluid height h and the local streamwise flow rate q. In the Newtonian case, a 
first model, also called integral boundary layer model was derived by Shkadov fiS] , by 
assuming that the flow remains close to the Poiseuille flow. However this model is not 
consistent at order one with respect to e. This results in an error on the prediction 
of the critical Reynolds number above which constant states are unstable under large 



tions, an information that is importan^to^p 
:ently, Ruyer-Quil and Manneville IHtpf R 



redict the onset of 
RQ2] corrected this 



wavelength perturbat 
roll-waves. More rece: 
shallow water model by using a Galerkin type approach and even derived consistent 
second order models by introducing a third variable r that measures the departure of 
the wall shear from the shear predicted by the semi parabolic velocity profile (which 
is precisely the quantity that has to be determined in order to find consistent models). 

There are many attempts found in the literature to derive a shallow water model 
for a thin layer of power-law AuM d^vn^ ^jj^^g^ plane, all based on the Shkadov 
approach and thus inconsistent l[NM, B*HS,' DM1, DM2, HCWL]. This is all the more 
problematic that these models are used to study the linear stability of power-law 
films since a direct Orr-Sommerfeld analysis is not available. Indeed, as pointed out 
in IRQ], the divergence of the apparent viscosity prevented (especially in the case 
< n < 1), up to now (see thereafter), to write a suitable linearization of Navier- 
Stokes equations about a Nusselt type solution. A first order consistent shallow water 
model was derived by Amaouche et al l[ADB], using the method of Ruyer Quil and 
Manneville. The first order models obtained in this way are not unique and moreover 
are not conservative, which can cause trouble^n the presence of shocks. In the mean- 
time, Fernandez- Nieto and the two authors l[FNV] proposed an alternative model that 
is first order and conservative by a more direct approach based on the expansion of 
the fluid velocity ^nd pressure field up to order one with respect to e. Very recently, 
Ruyer Quil et al |RQ] proposed a second order shallow water model in order to take 
into account of stream wise diffusion of momentum that plays an important role in 
the interaction of solitary waves when steady flow is unstable: here a regularization 
procedure was needed to apply Galerkin type method when n < 1 and various incon- 
sistency appeared in some terms found in the model. 



The purpose of this paper is two fold. On the one hand, we aim at deriving consistent 



1 INTRODUCTION 



4 



shallow water models that are consistent up to order two. On the other hand, we want 
to validate rigorously this computation by comparing the spectral stability of con- 
stant states of the shallow water model and the spectral stability of the corresponding 
constant states (Nusselt flows) solution of the Navier -Stokes equations. Note that 
up to now, there was no such direct study of the spectral stability of constant states 
for power law fluids down an incline since Orr Sommerfeld equations could not be de- 
rived. The main issue here, as pointed out in IKQ], is to take care of the divergence of 
the apparent viscosity. Indeed, if n < 1, one cannot linearize r = fi p \\D(u) || n-1 D(u) 
about the Nusselt type solution when D(u) = 0, i.e. at the free surface. Since the 
problem is only at the free surface and for the apparent viscosity, we have chosen 
to keep fluid strain rate as a (bounded or integrable) variable. Moreover, we intro- 
duce a (new) weaker form of the Cauchy momentum equations by integrating the 
Cauchy momentum equations with respect to the cross stream variable on an interval 
(z,h(x,t)) for all (x,t). In that setting, we require less regularity for the solutions 
in the cross stream variable. Indeed, instead of assuming that solutions are bounded 
with bounded derivatives, we only need that solutions lie in L p ((0, /i),R m ) functional 
spaces defined as follow 

rh(x,t) 

L p z ((0,h);R m ) = {v : Q t ^ R m | V(x,t) / \v(t, x, z)\ p dz < oo}, 

Jo 

where Q t denotes the fluid domain at time t. As a consequence, our new formulation 
does not contain the free surface boundary conditions and we can formulate a gen- 
eralization of Orr Sommerfeld equations to non Newtonian power law fluids n ^ 1. 
Moreover, we can propose with this approach a unified treatment of the shear thinning 
case n < 1 and shear thickening case n > 1. Though, depending on the case treated, 
we have chosen the classical constitutive law if n > 1 and the inverse constitutive law 
if n < 1. We have carried out a computation of the dispersion relation in the large 
wavelength regime up to order 3 which is important to decide whether second order 
shallow water models are consistent: in particular, we recover and quantify the fact 
that the model derived in fRQ] is not consistent except for n = 1. 

In order to derive consistent second order shallow water models, one has to com- 
pute an expansion of the velocity field up to order two with respect to e. This expan- 
sion is valid in a L p sense only under the assumption that n > 1/2. This threshold 
for n is^portant since it seems related to conditions for Navier-Stokes to be well 
posed IMNR]: in particular, for some applications in ice sheets dynamic (n = 1/3), 
one has to introduce a regularization of the power law and the linear stability and 
shallow water models that may result should depend drastically of the regularization 
process. Note that, in contrast to IRQ], in our approach, no distinction between the 
shear thinning case (n < 1) and the shear thickening case (n > 1) is needed. 
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I sec2 

The paper is organized as follows: in section 12, we introduce the set of equations 
that describes the motion of a thin power-law film down an incline aruigWrite it in 
a non dimensional form and in the shallow water scaling. In section |3, we move to 
the derivation of shallow water equations: we first compute an expansion of the fluid 
strain rate, pressure, velocity up to order 2 and then deduce a second order shallow 
water model. Throughout the study, we have considered the natural no-slip condition 
at the bottom. Other boundary conditions can be considered as well; like the Navier 
precondition which is usee g to 4 study thin film flow down a porous inclined plane 
UMBR] . Finally, in section k, we derive the Orr-Sommerfeld equations and compute 
an expansion of the dispersion law in the small wavenumber regime up to order 4. 
We then compare with results obtained for various shallow water models. Finally, we 
conclude and address some perspectives to this work. 



se 



^1 2 Governing equations 
2.1 Power fluid constitutive law 

We find in the literature various parametrizations of power law fluids, which depend 
on the possible applications. Usually, one describes a power law fluid (denoted PL in 
the sequel) as a generalized Newtonian fluid. Then the constitutive law reads 

^T^l (2-1) r = /i p || J D(u)|r t - 1 J D(u), 

where u = (u,w) T is the fluid velocity, D(u) is the deformation tensor and r is the 
fluid strain rate and D(u) is the deformation tensor which is defined as 



D(u) 



2d x u d z u + d x w 
d z u + d x w 2d z w 



|2 _ rr / p./ \ n /..\T 



whereas ||D(u)|| = Tr(D(u)D(u) )/2. This law is known as the Ostwald-de Waele 
or Norton Hoff and was first formulated for applications in metallurgy. Note that 
in the shear thinning case n < 1, the apparent viscosity /x p ||£)(u) || n_1 diverges and 
becomes unbounded at zero deformation rate. In order to deal with this singularity, an 
alternative approach consists in inverting the constitutive law so that the deformation 
tensor is a function of the fluid strain rate. Then, one Imp JJr|| = // p ||D(u)|| n and 
||D(u)|| = (yUp) -1 /™!^!! 1 /™. By inserting this relation into (BTTjJ one finds 



pl^TI (2.2) D(u) = (^-^Hrll^-V. 



One recovers the Glen-Nye formulation of the constitutive law for PL fluids by setting 
A = (/ip) -1 /™ and n* = 1/n: 



gl (2.3) D(u) = AWtW^-W. 
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Usually, the coefficient lies in (2, 4) and fdj most glaciers n* = 3 which corresponds 
to n — 1/3. We use the parameterization (tl.ff of PL fluids to write the Cauchy mo- 
mentum equations whereas we use (a rescaled version of) the parameterization (Llf) 
to carry out asymptotic expansions and to linearize Cauchy momentum equations. 



nsud 
nswd 
nsincd 



bcld 



bc2d 



impd 



2.2 Navier- Stokes eqns for a thin PL film down an incline 

We consider an incompressible power law fluid, with constant density p, flowing down 
an inclined plane, with a slope 9, under the action of the gravity (here g denotes the 
gravity acceleration). The two dimensional Navier-Stokes equations read, in the fluid 
domain Q(t) = {(x, z) | < z < h(x, t)}: 

(2.4) p (d t u + ud x u + wd z u) = pg sin(0) + d x (r xx - p) + d z r xz , 

(2.5) p (d t w + ud x w + wd z w) = -pg cos(0) + d x r xz + d z (r zz - p), 

(2.6) d x u + d z w = 0. 

This system is completed by the following boundary conditions at the free surface 
z = h(x, t): 

d xx hd x h 



(2.7) 



-d x h{j xx —p)\ x =h(x,t) + Txz\z=h(x,t) = 



l + d x h 2 fl 21 



(2.8) 
(2.9) 

whereas at the bottom z 



d x hr xz \ z=h ( x j) + (t zz -p)\ z =h(x,t) 



= a- 



d xx h 



(1 + d x h 2 ) 3 / 2 ' 
d t h + u\ z=h ( x j)d x h = w\ z=h ( x j), 

0: u\ z=Q = w\ z=0 = 0. Note that we take account of the 
surface tension: the capillary coefficient's denoted a. The fluid strain rate r is related 



surlace tension: trie capillary coemcie^nt^s denoted o. 1 ne nuid strain rate r is related 
to the deformation rate either by (Iz Tfor 0? Let us now write this system into a 
non dimensional form and in the shallow water scaling. For that purpose, we assume 
that the characteristic depth of the film is H, U its characteristic velocity and L is a 
characteristic longitudinal wavelength. The shallowness of the film is determined by 
the aspect ratio e: in what follows, we will assume e 1. The variables (x,z,t) and 
unknowns (r,p,u,w,h) scale as: 



x 



Lx', 



Hz', 



L / 



U, 



u = Uu', w = eUw', p = pgHp', r = fj, p (— )V, 
whereas the deformation rate scales as D(u) = ^D'(u') with 



/ 2ed x u ()■ ii 
V{U} - [d z u + e 2 d x w 



e 2 d x w 



2ed z w 
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We further introduce, respectively, the non dimensional Froude number, the Reynolds 
number and the Weber number: 



El 



Re 



P H n U 



njjl—n 



w e . 



pHU 2 



Then, by dropping the ' and s £tt^Dg 2 A =.i2 e sin(6')/F 2 , one may write 
with the boundary conditions (b.V,E.8,l.y) as 



Insufliswdh. 



sincd 



1 nsua | (2.10) 



d z (r z 



ReP 
F 2 



d z r xz + A = eR e (d t u + d x u 2 + d z (uw)) + ed x ^ 
) — A cot(#) = —ed x T xz + e 2 R e (d t w + d x (uw) + d z w 2 ) , 



ReP 
F 2 



d z w + d x u = 0. 

tolgni 

The power law (12.2) reads 



1 — n 



l-n 



Plgnia| (2.11) 2ed x u = (r 2 x + r 2 z ) 2n TxX) d z u + e 2 d x w = (r 2 x + r 2 z ) 2n r xz . 
whereas = — r rr . 



bcla 



Next, we set k = e F /W e : in order to see the influence of surface tension at 
%cid^c5^m cptli respect to e, we assume that k = 0(1). The boundary conditions 
(b.V,5.8,E.y) at the free surface are given by 

/ 010 v | / ReP,, o i eR e d x hd xx h 

(2.12) r xz \ z=h - e[T xx - -j^)\ z =h(x,t)Oxn = 



ReP, 



F 2 {l + (ed x h) 2 f/ 2 ' 



bc2a l ( 2 - 13 ) (Tzz- -j^)\z=h(x,t)-£T xz \ z=h(X:t )d x h = K— 2 / 1 + yn 2 )3/2' 

bc3a| (2.14) d t h + u\ z=h ( X! t)d x h = w\ z=h ( Xit ), 

whereas the no-slip condition at the bottom is unchanged: u\ z=0 = w\ z=0 = 0. 

Now we search for an Hilbert expansion of the fluid strain rate, velocity field with 
respect to the aspect ratio e < 1: since the apparent viscosity diverges at the free 
surface, this expansion does not make sense in Lf(Q,h(x,t)). Though this singular- 
ity is integrable and expansions will make sense in Lf(0, h(x, t)) spaces under some 

restrictions on n G (0, 1) and for suitable choices, of n n£ [1, oof. Hence we introduce a 

giswa I n sine a bcla pc2a L fbc3a L 

weaker formulation of (12.10, 12.11, 12.11) and (12.12, 12.13, 12.14) that is adapted to this 

l ns R^ pswa 

new functional spaces. This weak formulation is obtained by integrating (ETUI, E27TT) 
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on the interval (z,h(x,z)) and by using the boundary conditions (i.12, b. 13) . One 
finds 



r xz = X(h - z) - 



KSRe 



d x hd xx h 



+ ed x 



nsuf 



(2.15) 



uw 



+ 



ReP 
F 2 



F 2 (1 + {ed x h) 2 f/ 2 
sR e ^d t ( J udz) + d x ( J u 2 dz) 

" ' ' ' '' e3 ( f 

J z 



R ePi~ 

xx p 2 az 



:\cot(0)(h-z) - 



t xz dz) 



nswf | (2.16) 



+e 2 R e d t 



F 2 (1 + (ed x h) 2 fl 2 

rh r-h 

I w dz) + d x ( I uw dz) — 

J z J z 

i nsug. nswa I bcla [be 2 a , 

Any smooth solution of (12. iu, 12. ii) and (12.12, 12.13) is obviously a solution of 
(E7fe,E7T6"). Moreover, this weak formulation still holds for (r xx , T xz ,p) G L\{QJi^x, t)) 
and (u,w) G L 2 Z ((0 , h(x , t)) . Note that one can also write a weak form of (B7l?) by 
integrating the divergence free condition on (z, h(x,t)): one finds 



cmassz | (2.17) 



w = d t h + d x 



u dz). 



One purpose of this paper is to derive shallow water equations for h and q — f^u 
from the Navier-Stokes system in the shallow water scaling. We find an exact sys- 

tasuf i cmassz 

tern of evolution for [h, q) by setting z = in (12715) and (12717) and by using the 
impermeability condition w \ z= q = 0: 



Cmass | (2.18) 



d t h + 9 X 



%( j udz) = 0, 
d t (J udz) + d x {^J u 2 + j^j = (A/i - t xz \ z=0 ) + d x (j^ 



T~ X c 

~R f 



dz) 



Cmomentum (2.19) 

The mass cons 



K 



d x hd xx h 



F 2 (1 + (ed x h) 2 fl 2 ' 

Cmass 

rvation^law (12.18) is written in a closed form with (h, q) whereas one 
needs to write (t2.19 J in a closed form. This latter step is done by 

1. Computing an asymptotic expansion of r xz ,r xx ,u,w,p with respect to e, 

2. Computing an expansion of q = u dz with respect to e, 

3. Writing J Q h u 2 , p, r xx and t xz \ z=0 as function of h, q up to a given order 
with respect to e. 

We will see that there are various choices for this closure step but some of them are 
more natural and various constraints may also reduce the possibilities. 
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3 Derivation of shallow water models 

In this section, we compute an expansion of the fluid strain rate, fluid pressure and the 
velocity field with respect to e. Then we deduce lubrication models that are useful 
to study the stability of constant states of the full Cauchy momentum equations. 
Finally, we derive a consistent second order shallow water model for power law fluids. 



3.1 Expansion of solutions to Cauchy momentum equations 
in the shallow water regime 

hsuf faswf bgindraLgnia 

We assume that (u,w), (r,p) solution of (ET15',E7TC,ETTT^27lT) expand 



as 



(u,w) = (u^\w^)+e(u^,w^)+e 2 (u^,w^) + 0(e 3 ), 
(r,p) = (T^,p^)+e(T^,pW)+e 2 (T^,p^) + 0(e 3 ). 

By collecting O(e ) terms respectively in 'one finds 

D „(0) K D 

7-W = X(h(x,t) - z), rff + = Xcot(e)(h(x,t) -£)- -^d xx h, 

t£> = 0, d z u<® = X 1/n (h(x, t) - zf/ n . 

As a consequence, the fluid pressure is hydrostatic: 

= cos(#) (h(x,t) — z) — nd xx h. 

The fluid velocity (u^°\w^) vanishes at the bottom and is given by 

nA 1 ^ 



«(°) = (h(x, t) 1+l ' n - (h(x, t) - z) l+l ' n ) , 

w^ = - J" d x u^dz 



= -A 1/n ( fe (M) 1+1/ra + -^—(h(x, t) - z) 1+1 ' n - h(x, tf' n (h(x, t) - z)) d x h. 
\ n + 1 n + 1 / 



As quoted in |RQ], the apparent viscosity (c^-u^)"" 1 diverges at the free surface; 
though this singularity is integrable and asymptotic expansions hold in our mathe- 
matical framework. 
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Note that one easily recovers the constant, Nusselt type, solutions of the Cauchy 
momentum equations by setting h(x, t) = h: 

u= 7 ^(h 1+1 ' n -(hh-z) 1+1 / n ), 
f xz = X(h - z), = A cot(9)(h - z), 



I XX 1 zz 



0, W = 0. 



Hence, by restricting the analysis to small amplitude motions \h(x,t) —h\<^l, we 
expect that the asymptotic expansion shall provide a useful information on the stabil- 
ity of constant states under large scale perturbations: this is done in the last section 
of the paper. 

Next, we collect 0(e) terms in (6fS,6?S) and (t^fljf 





xz 


= ~Re 


(3.2) 




-d x ( 


(3.3) 
(3.4) 


ReP (1) 

F 2 
r (i) 

XX 


= -d x ( 

V xz 



u (0) w (0) 



I e 2 dz) - —^d x hd xx h, 

J Z 



T (o)\ _ T (i) 

XZ I XX 1 



t<3 =2(TW) 1 - 1 / n d x iiW = 2A (h 1/n (h - zf- x l n -(h- z)) d x h. 



Hence, the correction to the fluid pressure is given by 

= A ((h - z) - 2h l ' n (h - zf- l ' n ) d x h. 



ReP (1) 



F 2 

Insuf |nswf 

The weak formulation (Em.ErnT) holds only if (r,p) G ^((0, h); M 2 [ Jt J ; x K): wc 



easily deduce from the expression of t^x that our calculations only make sense if 
n > 1/2. In what follows, we will see that there is no additional restriction on n in 
order to compute an expansion of the fluid strain rate and velocity field up to second 
order with respect to e. Note that such a restriction on n is also met when proving 
the well posedness of the Cauchy problem for power law fluids in bounded domains 
with periodic or Dirichlet (no-slip) boundary conditions IMNR]. 

Now we can compute t^z (and thus u^). Note that it contains time derivatives in h: 
it is standard in the derivation of thin film models that one can exchange time and 
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spatial derivatives. More precisely, by using the mass conservation law (12.18) which 
reads 

rh 

d t h + d x ( u {Q) dz) = 0{e). 
Jo 

and by using the expression for u^°\ one finds 

so that one can replace dth by —^^h^^^dji in the expression for t^z ■ A lengthly 
but straightforward computation that was implemented in a formal computation soft- 
ware provides the following expression for t X z : 

T a, = (2» + -*) - n(h - 

(2ra + l)(ra + l) 

+n-^(h - z)d xxx h - Xcot(9)(h - z)d x h. 
In particular, the tangential bottom strain rate expands as: 

r^| 2=0 = ri° } |, =0 + st^\ z=q + 0(e 2 ), 
1 J r^| 2=0 = \h + e I - 2n+1 ~ A/i oot(fl) J 9 x /i + 0(e 2 ). 

We complete the description of the velocity field at order one by computing succes- 
sively dzU^Kuty = J Q Z dzii^dz and = — d x u^dz. Indeed, by collecting 0(e) 
term into ( , one finds 

d z u^ = U\{h-z)f^r^ 

and u^\w^ easily follow, which completes the description of the velocity field up 
to order one. For the sake of conciseness of the presentation, we only provide the 
expression of interest to us, namely the correction to the discharge rate q^: 

r h 

q W= / (h- z)d z u {1) dz 



( 2\^h^R e \^h^cot(6) \ nR e f 2+1/ 

V(2n + l)(3n + 2) 2n + 1 J ° x + (2n + 1)F 2 A ° xxx ^ 



so that the full discharge rate q expands as 
nh 2 (\h) l ' n 



2n+l 



+ e 



f 2X 3/n ~ 1 h 3+3/n R e _ A 1 /"/i 2+1 /"cot(^) \ 
V(2n + l)(3n + 2) 2n + 1 J x 



+ e , KR \ ^ \ l/n - l h 2+1 ' n d xxx h + 0(e 2 ). 
(2n+l)F 2 v ; 
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Note that t^J does not play any role in the first order expansion of both and 
Tx^\ z =o and thus there is no restriction on n G]0, 1[ to write^rsi order models. This 
is consistent with a previous paper by two of the authors FNV]: the expansion was 
carried out in Lf(0, h) which resulted in a complicated treatment of expansion near 
the free surface where the apparent viscosity diverges. Here, our weaker L\(Q,h) 
formulation of the problem provides a simpler approach to compute a first order ex- 
pansion of the solution and derive first order shallow water models. 

In contrast to (fNV] , we are now in a position to derive second order expansion 
Snsu? rliwi^i^ilP^ 3 ^ s ^ ram ra ^ e an< ^ discharge rate by using the new formulation 
of the Navier-Stokes equations under the restriction that n > 1/2. 
Collecting 0(e 2 ) terms in (&.fl), one finds 

t£> = 2(X(h - z)) X - Xln d x vF> + ^(X(h - z))-*t<£ r« 

Note that this function is Ll(0,h) only if n > 1/2 but t^x does not play a role in 

(2) 

the derivation of a second order shallow water model. Next, we compute TxJ which, 
in turn, provides an expression for and q^ 2 K One obtains t£) by collecting 0(e 2 ) 
terms in (E7T5): 

rf} = -R e (dt £ + 2d x £ A« - (u^wW + u^w^ 

0.6) + ^\«_^" 

The first part of the expression of ri^ which lies in front of R e is easily proved to be 
L™(0,h) whereas the last part is O ((h — z) 1 ^ 1 ^)) which is L\{Q,h) provided that 
n > 1/2. Next, by setting z = in the expression of , we obtain the second order 
correction to the bottom stress: 
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2A^~ 1 /^+ 4 (5n 2 + Un + 6) 



1 (2n + l) 2 (3n + 2) (4n + 3) 1 



com i f Al ^: 3 / n+2 Q <w> + xh ] i2l '- :]) »,j, 

1 (2n+ 1) (3n + 2) / 



+Ri 




, , 2 A- (n + 1) (12 n 2 + 46n + 21) , :i , 

n(2n+ l) 2 (3n + 2) (An + 3) 

(n + l)^ 2 /"" 1 ^/™ 2 /(l-n)^) 4 6h(d x h) 2 d xx h 3h 2 (d xx h) : 

-ti p 5 1 7, V 



(2n+l)F 2 e V n 3 n 2 n 

reA 2 /"-^ d2 / (2n + 3)h 2 d x hd xxx h (n + 2)h 3 d xxxx h 
+^ ™ -^C 7^ — — ^ H 



F 2 e V n(2n + l) (2n + l)(3n + 2) / ' 

Now, we compute successively J^4^?\ u ^ — Jq d z u^dz and = J^(h — z)d z u^dz. 
By collecting 0(e 2 ) terms in (ts.fi), one finds 

d z uV = - dxW l°) + I (r (0)^ " \g) _ ^(r£>)^ " 2 + (r«) 2> 

A careful inspection of the various terms in the expansion of d z u^ shows that 
\d z u^\ = 0((h-z)~ l / n ) which is not L*((0,/i)). However, = f*d z u& = 
0((h — z) 1 ^ 1 ^ n ) is integrable on (0, h) and the second order correction to the dis- 
charge rate q defined as q^ = J^u^dz is finite. A lengthly but straightforward 
computation yields the following expression for q^: 
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= / 4Ai-'^( 8n ' + 2 5n+ 12) \ 
e \{2n+lf (3n + 2) {An + 3) {5n + 4)J 

- RcCom ( 2A*-ft*"(5n» + 14n + 6 A 
V ; \(2n + l) 2 (3n + 2) (4n + 3)y 

4A^+ 3 (n + 2) , 

H ; ~~t o xx h 

3 (2n- 1) (3n + l) 

2 f 3 A"~ 2 fri+ 4 (n + 1) (35 n 2 + 226 n + 120) \ 2 
e \ v 2n(2n+l) 2 (3n + 2) (4n + 3) (5n + 4) ) x 

, / r (d x hf 

n (2ra + l) (3n + 2) y v ; 

A" /^+ 2 (10n 2 + 25n + 7 - 3(n - l)(2w - l)cot 2 (fl)) 2 
6n(2n-l) (2n+l) ( x } ' 

K 2 Rl{l-n){\hf/ n h\d xxx hf , KR e (Xhy/ n (n-l)h 2 d x hd xxx h 



F A 2n(2n+l)A 2 F 2 n(2n+l)A 



cot(0) 



| K(n+l)Rl(\hfl n / 2(1 -n)h(d x h) 4 + I2h\d x hfd xx h + e/r 5 ^/^ 



(2n+l)A 2 F 2 V (3n + 2)n 3 (3n + 2)n 2 (3n + 2)n 



(2ra + 1) A 2 F 2 V n(3n + 2) (2n + 1) (3n + 2) (4n + 3) 



This completes the description of the fluid strain rate and velocity field up to order 2 
with respect to e. We will use these expansions to provide various thin layer models 
both one equation (lubrication) models and two equations (shallow water) models. 



3.2 One equation-lubrication models 

In the previous section, we expanded the fluid strain rate and velocity field with 
respect to the aspect ratio e in the shallow water scaling e — > 0, where the coefficients 
are functions of the fluid height h and its spatial derivatives and the cross stream 
variable z. As a consequence, these quantities are completely determined if the fluid 
height h is known. One derives an equation for the fluid height by considering the 
mass conservation law 



d t h + d x q = 0, q = 
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and by substituting the expansion of q with respect to e. Indeed, by inserting an 
expansion of q to order j = 0, 1, 2: 

jubn] (3.7) d t h + d x q^ = 0(e j+1 ), j = 0,1, 2, 

and dropping 0(V +1 ) terms, one obtains a hierarchy of lubrication models consistent 
with the full Cauchy momentum equations and of order j. The first model of this 
hierarchy is a Burgers type equation 

l^bOl (3.8) d t h + d x - 0, 

which exhibits singularities in finite time. Though, it provides a useful information, 
namely the phase velocity of small perturbations of constant states c(h) = (\h) l l n h. 
The next model in the hierarchy, obtained for j = 1 is a generalization of the Benney 
equation: 

/ nh 2 (\hy/ n \ _ // A 1 /"/i 2+1 /"cotan(g) 2A 3 /"- 1 fr 3+3 /"i? e \ 
* + * V 2n + l )~ £ x \\ 2n + l (2n + l)(3n + 2)) x 

D \l/n—l 

( 3 - 9 ) -^ (2 J +1)F2 ^(fe 2+1/n ^^)- 

l lubl MM IFNV 

Equation (t5.9) was already derived formally in |NM] and in l[FNV] . In the former 
paper, the au^hop proposed two equations shallow water models that are na^on- 
sistent wittL^g^T (even in the Newtonian case). As already noted in fNM], |FNV], 
the model (l3T9j yields a stability criterion of constant states under small wavenumber 



Re_c 



perturbations. More precisely, by linearizing (l3.y) about the constant solution h — 1, 
one finds 

dth + }}l n d x h = — d x ( (cotaxi{6) — — —)d x h] — ke-^— — r-=^d xxxx h. 

2n + 1 V 3n + 2 J ) (2n + l)F 2 

Then, the constant state h — 1 is stable if the following criterion is satisfied 

3n + 2 

(3.10) i? e < ^^cotan^) = .Ree. 

If the meaji velocity of the Poiseuille flow is chosen as the characteristic velocity, 
equation (BTO) is a generalization of the classical stability criterion for Newtonian 
fluids: 

5 

Re < -cotan(6 ) ). 



3 DERIVATION OF SHALLOW WATER MODELS 



16 



cmom n 



[Re c l lubl 

If (&TTTT) is satisfied, the equation (K9J~ is globally well posed for small amplitude 
solutions ^dii^J 1 decrease asymptotically to 0. However, if this criterion is not satisfied, 
equation (1379) is ill-posed and one has to consider a higher order model. The next 
model in the hierarchy is 

dth + d x (qV>(h,d x h,d xx h)) =0. 

I lubl 

However, the next correction term here is dispersive and would not stabilize (13.9): one 
should then consider a third order model that would yield a Kuramoto-Sivashinsky 
type equation. However, this type of equations only makes sense near the instability 
threshold < Re — Re r <C 1 which is limited for applications. ,, , 

Jlubn 

In the last section of this paper, we relate the lubrication models (1377) with the 
linear stability of constant states of the full Navier-Stokes equations by a direct anal- 
ysis of the linearized Navier-Stokes equations, so called Orr-Sommerfeld equations 
(which are written here in a non classical form). 

3.3 Two equations- Shallow Water type models 

In this section, we propose an alternative to lubrication models and introduce two 
equations shallow water models with a source term that 

l lutaO I lubl 

1. contain the lubrication models (BTSj and tnHn the limit e — > if Re < Re c , 

2. are locally well posed if Re > Re c for order one models, 

3. are globally well posed for second order models. 

,„ The two equations-shallow water models are derived from the integral equations 

gmagCmpffgntum Cmasg 

(12.18,12.19). Equation JI^TS) is exact and already written in a closed form. There 
remains to write (LlyJ m a closed form: 

d t q + d x ([ u 2 + ^-dz) = ^—(\h - t xz \ z=0 ) + d x ( [ r xx dz) - 



o 



F 2 eR e y XZIZ=UJ ' XK J xx ' F 2 (1 + (ed x h) 2 ) 3 / 2 



closl 



One obtains a order j + 1 model by considering the truncated equation 
(3.11) 

with 

j 

(3.12) (u ij \f (j \p ij) ) = ^e\u {j \T {j \p^). 

k=0 
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Then these functions are written as function of (h, q) and their spatial derivatives: in 
this process, various choices are possible. For the inertial term, one writes 

H^2\ (3.13) / (u^fdz = (3{n)\ + l(h,e,n) + 0(e j+1 ), 

Jo h 

where the form of / depends drastically on the choice of (3(n). The fluid pressure 
being hydrostatic, one has 



(3.14) I p j {z)dz = P{h,n,e) + 0{e j+1 ). 

Jo 

MM 

Finally, we write the bottom stress as it is found in the literature e.g. in [NMJ: 

EUsH (3.15) ?£ + %=o = ( (2 ^ 2 1)g ) n + *XM,e) + 0(e^ 2 ). 

According to (fNV] , one can write F (h,n,e) = d x f(h,n) for the first order model 
(j = 0). Moreover if we add the Galilean invariance as a constraint, one has to choose 
/3(n) = 1. This imposes the following form of the first order model: 

d t h + d x q = 0, 

(<f_ n 2 X 2 / n h 3+2 / n (3n + 2)Acot(fl)h 2 > 
'''' ' " ' ' 1 J + (3n + 2)(2n + l) 2 + 4(2n + l)R e 

3n+2 f xh-( q 2 -^ r 



2e(2n+l)R e \ ^ nh 2 

(3.16) +4 {3n + 2 \ hd xxx h. 

v ) F 2 2(2n + 1) 

This model is, up to 0(e 2 ) terms, equivalent to the first order models found in I(aDB] 
and IKQ] (in this latter paper, this is obtained by dropping streamwise diffusion 
terms) : 

Un + 6 qd x q 2n + lg 2 (3n + 2) 

9tq + -4^— ~ 3 ^T3h^ h + 2WTWe XC ° mhdxk 

3n + 2 / ,(2n + l)q, 
1 Ah — (- 



2e(2n + l)R e \ v nh 2 

k (3n + 2) 

sw_nc I (3.17) H - — -ho xxx h. 

1 v ) F 2 2(2n + 1) 
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Note that this latter system is not written in a conservative form. Indeed by using 
the expansion of q: 

q = (A/*) 1 ^ 2 + 0(e) = f(h) + 0(e) 

in the convection and source terms of and up to order 0(eR e ), one can derive 

a general family of shallow water models in the form 

d t h + d x q = 0, 

d t q + a(h, q)d x h + b(h, q)d x q - k(h, q)d xxx h = —— (g(h) - ip(h, q)) . 

eri e 

that are all consistent up to order one if the functions a, b, k, g, ip satisfy some con- 
straints at q = f(h). We will write these constraints explicitely in the derivation of 
shallow water models which are consistent up to order two. 

Up to now, there was only one attempt to derive a second order shallow water 
model for thin power law film [RQ]: however, on the one hand the model derived is 
not consistent with the case n — 1 for fully nonlinear term whereas we shall prove 
in the next section that quasilinear ones are also not consistent except in the case 
n — 1: indeed, as already mentioned in |RQ], this is due to the fact that correctors 
in the inertial terms are neglected. 

We set f3 = eR e and search the second order shallow water model in the form 

d t h + d x q = 0, 

d t q + a(h, q)d x h + b(h, q)d x q = — (g(h) - ip(h, q)) + nk(h, q)d xxx h 

+ P (P^h, q)d xx h + P 2 (h, q)(d x hf + P 3 (h, q)d xx q) 
+ K0 {N 3 , 2 (h, q)(d xxx h) 2 + N 1A (h, q)(d x h) A ) 
+ K0 {N 12 ,2i(h : q)d xx h(d x h) 2 + N 2j2 (h, q)(d xx hf) 
+ k(3 (N A) i(h,q)d xxxx h + N 31M (h,q)d xxx hd x h) . 

In the limit — > 0, this model provide an expansion of q in the form 

The shallow water model is consistent with the Cauchy momentum equations up to 
order two with respect to e if the following constraints are satisfied: 



g(°)=g(°), q U=R e qW, q <®=Rl<p\ 
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The first constraint = implies g(h) = ip(h, f(h)) with 

f(h)=q^\h) = ^(Xh)^h 2 . 

Next, the constraint q^ = R e q^ is satisfied only if: 

f{hf - a(h, f(h)) - b(h, f(h))f'(h) (Xh) 3/n h 3 {Xh)^h 2 



ip q (h,f(h)) (2n + l)(3n + 2)A R e (2n + 1) 

k(h,f(h)) (\hf/ n h 



cot(0), 



*l> q (h,f(h)) (2n+l)F 2 A' 

Finally, we consider the constraint q^ = R 2 q^ and A{h) = ip q (h : f(h)). First, we 
identify the term in front of d xx h in q^ and cf 2 ^: 

Pi + (Xh)^ n hP 3 {2h{Xh) 1 / n - b)((Xh) 2 ' n h 2 - a - (Xhf^hb) 
A + A 2 

(8n 2 + 25n+ l2)h 5 {Xhfl n 



A 2 (5n + 4)(4n + 3)(3ra + 2){2n + l) 2 
(5n 2 + Un + 6)(Xhf/ n h 4 
R e {4n + 3)(3n + 2){2n + l) 2 A C ° t( ' 
4(n + 2)(Xh) l ^ n h 3 
+ 3(3n+l)(2n-l)R 2 e 

Next, by identifying the terms in front of (d x h) 2 , one obtains 

nft + (n+l)(AA)Vft + q + ^)V--tf(^ (gA2(Afc)2/ „ _ &) ^ 

3{n + l)(35n 2 + 226n + 120)/i 4 (A/i) 5 / n 



2A 2 (5n + 4)(4n + 3)(3n + 2)(2n + l) 2 
■ (A ' !)3, "' i3 -cot W 



An(3n + 2)(2n + l)i? e 

^cot 2 (#) 



(A/i) 1/n /i 2 (-6n 2 + 9n - 3) , 



Qn(2n + l)(2n-l)R 2 e 
(A/i) 1 /™/ l 2 (l0n 2 + 13n + 19) 



6n(2n + l)(2n- l)i? 2 
Then, by identifying the term in front of (h xx ) 2 in q^ 2 \(f 2 \ one finds 
N 2>2 1 / (n+ l)(Xhfl n h 3 (n + l){Xh) 2 / n h 2 

TTT~/ T~, r <^ 



F 2 V nA 2 (3n + 2)(2n + 1) n(2n + 1)AA 
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By collecting the term in front of (h x h xxx ) in q( 2 \(£ 2 \ one obtains 
#31,11 2 (Sn + 7)(\hf/ n h 3 | (n-l)(\hy/ n h 2 cot ^ 



A nA 2 (3n + 2)(2n+ 1)F 2 nXF 2 (2n + 1)R, 
+ 



n + 1- 2n 2 )a + (5n - 2n 2 )bh{Xhfl n + (2n 2 + n/i 3 & 9 - (9n + 7)/i 2 )(A/i) 2 / n 



n(2n + l)\AF 2 

Then, we collect the term in front of (h x ) 2 h xx in q( 2 \ q^>: 

iVi 2 , 2 i = J_ / (n+l)(Afr) 3 /"fr 2 _ 6 (n + l)(A/i) 2 /^ 



A F 2 V n 2 A 2 (3n + 2)(2n + 1) n 2 A(2n + l)AF 2 
Next, we identify (c^/i) 4 in front of g (2) ,g (2) : 

iV M 1 ( (1 - n 2 )fr(Afr) 3/ " (n 2 -l)(A/i) 2 / n 



A F 2 V n 3 A 2 (3n + 2)(2n+l) n 3 AA 
Finally, we identify respectively (h xxx ) 2 and d xxxx h terms in q^ 2 \q^\ 
#3,2 _ « ( (n - l)(\h) l l n h 2 h\\hf/ n ^ qq 



A F 4 \ 2nA 2 (2n+l) 2(2n + l) 2 AA 2 y ' 

and 

N Ayl 1 ( n (5n 2 + 14n + 6)h 4 {Xh) 3 / n [b - 2h{Xh) 1 / n )h 2 {Xh) 1 / n 



A F 2 V (4n + 3)(3n + 2)(2n+ 1) 2 A 2 (2n+l)AA 

We search for a shallow water model that is consistent and is in the simpler and 
classical form: for that purpose, we try to cancel the fully nonlinear capillary terms 
Nij which contain four derivatives. First, if one sets N 2t 2 = 0, then A = ip q (h, f(h)) 
is given by 



1 , 1 



x n-i, i+n(3n + 2) , /, nh 2+ nX 
A = X~h ~- — - — - = ip q \h, 



2 Tq \ ' 2n + 1 

One possible choice for ip q (h,q) is then 



2n (3n + 2) 



^ (h, q) = q n ~ L (2 n + 1)" n~ n+v \i 

and by integrating with respect to q, one finds the natural form of the friction term 
for power law films 
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The capillary coefficient k is written as 

(3n + 2)h 



k(h,q) 



2F 2 (2n+ 1) 



As a byproduct of our particular choice, we also get = Ai 2 ,2i = N± t 4 = 0. In 
order to cancel A 4j i, one needs 

nh 2+ «\n i i+n (lln + 6) 

o(/i, ) = nX^h « 



2n+l y 8n 2 + 10n + 3 

and consequently, a suitable choice is 

v lH! h 4n + 3 
The consistency of the shallow water model at first order yields: 

nh 2+ ^X^ 3X^h 2 ^n 2 12n 2 + 17n + 6 Xh 

a{ ' 2n + l } ~ ~8n 2 + 10n + 3 + 2(8n 2 + 10n + 3) P~" C ° t( } ' 

One possible choice is 

s 3(2n + l)g 2 3n + 2 Xh /n . 

We also get iV 31jll = 0. Finally, there remains to determine P*: one suitable choice is 

9n(n + l)(5n - 4)A- 1+4/n /i 3+4/ri 



4(2n + l) 2 (4n + 3)(5n + 4)(3n + 2) 



1 (10n 2 + 25n + 7)(3n + 2) (2n+l)g 
31 J P 2 12(2n-l)(2n+l)(n+l) 1 n/i 2 ' ' 
2n 2 (l - n)A- 1+4 /"/i 4+4 / n 



(2n + l) 2 (4n + 3)(5n + 4)(3n + 2) 



1 n(3n + 2)(9 + 10n-29n 2 -14n 3 )g (2n + l)g _ t 
l ' J P 2 12(2n-l)(2n+l) 2 (3n + l)(n + l)/i l n/i 2 J 
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Note that the first order consistent model found in pEDB] and [RQ] is precisely the 
one for which the correction terms due to surface tension vanish. The second order 
shallow water model reads 

lln + 6 ^ _ 2n + l g (3n + 2) 

tH An + 3 h An + 3h 2 2(2n + l)R e v ; 



2(2n+l)/?V / F 2 2(2n+1 

+/3 (P^d^h + P 2 (h)(d x h) 2 + P 3 (h)d xx q) . 



Note that Pi and P 3 , which play a role in the linear stability of constant states has 
2n — 1 in the denominator and only make sense if n > 1/2. Once this condition 
satisfied there is no distinction between the shear thinning case n < 1 and shear 
thickening case n > 1. 



4 Direct stability analysis 

In this section, we study rigorously the linear stability of constant Nusselt flows 
without any regular izat ion of the apparent viscosity. We show that the expansion of 
velocity field and fluid strain rate carried out previously in the shallow water /small 
wavenumber regime provide naturally an expansion of eigenfunctions and eigenvalues 
of the Orr-Sommerfeld equations in the small wave number regime. In particular, we 
show that the dispersion relation associated to linearized lubrication equations about 
the constant state h — 1 yields the expansion of the eigenvalue A ~ of the full 
problem with respect to the wavenumber k — > 0. This analysis is valid up to order 
2, the order of the discharge rate expansion carried out previously: for completeness, 
we have computed an expansion of the eigenvalue A ~ up to order four with respect 
to k in order to generalize the classical results for newtonian fluids. 



4.1 Rigorous formulation of the Orr-Sommerfeld equations 

In this section, we write rigorously the linearized Cauchy momentum equations about 

a Nusselt solution: in this case, the apparent viscosity diverges at the free surface. 

Though, by using our weak formulation and keeping the strain rate as an unknown, 

we remove this singularity. We start from the non dimensional Cauchv momentum 
insua tol gnia J jpclgWc2gfcc3a 

equations (l27TIT,l2.ii) completed by the boundary conditions (I27T2,I27T3',I27T3') and fix 
e — l 1 . For the sake of conciseness, we have neglected surface tension in this discus- 



1 At this stage, we do not focus on a particular wavenumber regime 
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sion. 

In what follows, we study the stability of the constant (in x and t) Nusselt solution 

h = 1, r 0<xz = A(l - z), r , xx = r 0>22 = 0, p = c(l - z), 

W = 0, C/ = ^ (1 - (1 - z) 1+1/n ) ■ 
n + 1 ' 

We linearize the Cauchy momentum equations about this constant solution 

Q x p 1 

(4.1) d t u + U d x u + UqW + = — [d x T xx + d z r xz ) , 

() z p 1 

nswl| (4.2) d t w + U d x w + y^ = — (d x r xz + d z r zz ) , 

(4.3) d x u + d z w = 0. 

The linearized boundary conditions at the free surface are written 

I \h —\ Tzz \ c ^ 

7~xz\z=l All, \z=l j-, \z=l 



nsul 



nsincl 



be lin 



p2lz=l R ^ z=i p2 , 



d t h H —d x h = w\ z=1 . 

n + 1 

The constant solution has a singular apparent viscosity: in order to write the lin- 
earized problem, we invert the relation r = \\D{u)\\ n ~ l D{u) as D{u) = ||r|| 1 / n_1 r. 
Then, we linearize this constitutive law: 

1 — n 1 — n 



pl-lin (4.5) (A(l - z)) n r xx = 2d x u, - (A(l - z)) n r xz = d z u + d x w. 



Since the velocity field is divergence free, there is a potential function so that 

u = —d z (p, w = d x cf). 

I nsuhswhsincl fol lin I bc lin 

The system (14.1,14.2,14.3) and (14.5) with the boundary conditions (14.4) are homoge- 
neous in time and space: we then search for solutions of this linearized system in the 
form 

(<j),a,p):=e^ k *(<P,a,p). 

Then, one finds 



ReP 

F 2 



r xx + Re(U' <p - (A + C/ O )0')) > 



R e p' 

~f2~ - T 'zz = ik ( T xz ~ lk R e {A + U )(/>) , 
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whereas the boundary conditions read 



tHI] ( 4 - 6 ) „\l/n 



'xzU=i = A/i, ^r|*=i - r^| 2= i = Acot(0)h, 



nA 1/ri 

x=1 = (A + ——)h, (p\ z=0 = <f>'\ z =o = 0. 
n + 1 

The linearized constitutive law reads 

1 — n 1 — n 



picll (4.7) (A(l-z)) n T xx = -2ik<l>', (A(l - z)) n t xz = <f>" + k 2 <f>. 

1 71 



Inscl 

Following the strategy of the previous section, we integrate (14.6) over the interval 
{z, 1) and use the boundary conditions (kib) : this yields 

nsulf | (4.8) T xz = \~h-ik^ (^-r xx + R e (U&-(A + U )(f)')^dz, 



nswlf | (4.9) ^0- + T xx = \cot(6)h-ik j r xz dz - k 2 j R e (A + U )<f>dz. 

" J z J z 



This comes with 
tjl] (4.10) r ra = -2^(A(l-^)) 1 - 1/ V, 



ihill (4.11) 0= - f f -(A(l - y)) 1 /"- 1 ^ + k?<f>(y)dydz. 
Jo ./o n 

In particular, one has 

<f>\ z =i = - f ^—^(1 - iQ 1/n T^ + A; 2 (l - 
Jo n 

4.2 Connection with the shallow water expansion 

In what follows, we will consider the low wavenumber regime A: — which corre- 
sponds to the shallow water regime for the Navier Stokes system: we shall prove 
that Orr Sommerfeld equations are a linearized version of the "weak" formulation of 
Navier-Stokes equations and that expansion of eigenvalues and eigenf unctions of Orr 
Sommerfeld are obtained directly by the linearized iterative scheme. 

In the shallow water scaling, r, u, w, p are functions of h and its partial derivatives 
and one can compute an expansion with respect to e up to order 2: 

2 

(r, u, w, p)(z, x,t) = Yl £k (r {k \u^,w^ k \p^)(z, h(x, t)) + 0(e 3 ). 

k=0 
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We denote 
(T ,Uo,Po) = (ri J,u^,p^)(z,l), 

T = A(l - z), U = n\ l l n {\ - (1 - z) 1+1 ' n )/{n + 1), ^ = Acot(0)(l - z). 

hsuf I nswf i cmasgz 

By differentiating the weak formulation of Cauchy momentum equations (E2.i5,Li6,LiY) 
and (I5.fi) about h — 1, one proves that 

(?,u,w,p) = d(r,u,w,p)\ h=1 (h), 

satisfies the problem: 

1 — n ] 1 — n 

2ed x u = T n f xx , d z ii + e 2 d x w = -T n f xz , u\ z=0 = 0. 

n 

and 

f xz = \h + ed x ^ jT f xa! - ||(p + P (l)^ 

-eR e (d t {j^ u + U (l)h) + d x (U 2 (l)h + 2 U Q u)-U Q w^J, 
f xx + ^- = \cot(0)h - ed x Qf f xz + T (l)~hj 
+e 2 R e (d t {j^ *>) + d x (£ W)), 

i cmasgz 

whereas the linearized mass conservation law (12717) is written as 
w = d t ~h + d x ^ u + U (l)1iJ , w| z=0 = 0. 
Then, it is a straightforward computation to prove that 

f x « = \h + ed x ( [ - ^ 



2 



-eR e (d t {J^ u) + ^(^ U u + ^ 
t X x + = \cot(6)h -ed x (J f xz ) 

+eR e (d t (£ ™) + dx(J z U w) 
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These equations are homogeneous in x and t: we search for solutions in the form of 
plane waves 

(f,«,TS,|0:=e iCAt+<Cx (r( 2! ),-^),^),p), 

which yields precisely the linearized Orr-Sommerfeld equations ( TTJ and 

(H~B) by setting /c = e£. As a result, one obtains directly an expansion with respect to 
k of eigenvalues (A, 0, r) simply by identifying the coefficient in the Taylor expansion 
with respect to e 

(4.12) A = A(°) + A;A« + k 2 A {i) + 0(k 3 ) 

(4.13) = A<°> + e£A« + e 2 eA {2) + 0(e 3 ). 

Indeed, the expansion of A is easily obtained through the mass conservation law 

d t h + d x {[ u + U o (l)h) = ^Ah+ [ u + U o (l)h = 0. 

J z JO 

Recall that q = u so that q — dq\ h= Ji — I u + U (l)h. Moreover, we have proved 

Jo 

in the previous section that q expands as 

q = q( \h)+eq^\h)d x h + e\q^\h)d xx h + r^(h)(d x h) 2 ) + O(e 3 ), 
so that q is given by 

~q = (qW)'(l)h + eq {1 \l)d x h + e 2 q {2 \l)d xx h + 0(e 2 ) 
From this expansion and by using the relation A + q = 0, one finds 



A<°> = A« = A< 2 > = <r 



,(2) 

21 — yd j w> 11 ~~ "H i 1S - ~ s 

This yields 



and 



A (o. = _ A V», A q, = _V VM»)\ 

V(2n + l)(3n + 2) 2n + 1 / 

A ( 2) = R 2 ( 4 A^ 2 (8n 2 + 25n + 12) 

(2n+ l) 2 (3n + 2) (4n + 3) (5n + 4) 

2\i~ 1 (5n 2 + 14n + 6) 



-R e cot(0) 



(2n+ l) 2 (3n + 2) (4n + 3) 



(4.14) 



4 A- (n + 2) 



3 (2n- 1) (3n+ 1) 



Similarly, one can obtain an expansion of the potential and fluid strain r up to 
order 2 with respect to e. 
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4.3 Higher order correction to Orr Sommerfeld 
eigenvalues/eigenfunctions 

We stopped the full non linear iterative scheme up to order 2 both from applications 
point of view (it is sufficient to derive second order shalllow water mdels) and due to 
some limitations due to the integrability of nonlinear terms in the expansion. How- 
ever, at the linear level, one can continue the expansions one step further. As an ap- 
plication, one could deduce a consistent amplitude equation of Kuramoto-Sivashinsky 
type near the instability threshold. 

Recall that we expand r, <fi,p with respect to k(— e£) as k — > 0: 

( r ,0,p) = ( r (°),0(°),/>)) + ^ 

Collecting O(k ) terms into (S S 8jl9) lf and (StO) yields 

R e p^ 



4°J=\h, T £>=0, = Xcot(9)h. 



We then deduce that 

= \ l/n ( — + (l-z)- — (1 - z) 1+1 / n ) . 

\ n + l n+l J 

This provides 

(o) U = A(°) + ^ = ^-A 1 /", 
1 n + l n+l 

so that A(°) = —A 1 /™. Next, ^e^ompute higher order corrections. For that purpose, 
we collect 0(k) terms into (k.») and (B.9): this yields 

T<£ = -i [ + Re(U^ - (A<°> + C/ O )(0 (O) )') d* 

= -i Acotan (0 1 - 2 + — — v - ; . 

V ra + 1 In + 1 n 7 

Then, one obtains 

<^-^(^-^^^-^>)- 

By integrating twice this equation with the boundary conditions 0| 2 =o = <f>'\ z =o — 0, 
one finds 

\2+l/n\ 



A 3 /™- 1 (n{l - (1 - z) 2+l l n ) n n(l - (1 - zf +2 ' n ) 



n + l) 2 V 2n + l 4n + 2 v 3n + 2 
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A ( } = 7 COt(0) — 

2n + 1 V 3n + 2 y 



This yields 

As a result, we obtain the following expansion 

\l/n / op \2/n— 1\ 

A(fc) = -ik\ 1/n - k 2 — cot(0) - =^ + 0(A; 3 ). 

w 2n + 1 V 3n + 2 y v y 

which provides the stability criterion under small wavenumber perturbations 

(3n + 2)A 1 - 2 /" 
i? e < cot(0). 

MM IFNV 

This criterion was already derived in |NM] and l[FNV] from a consistent lubrication 
theory. However, this is the first rigorous result that is derived directly from the 
Cauchy momentum equations through a Orr-Sommerfeld type analysis. The key point 
was to introduce the weak formulation so that one can obtain a well posed linearized 
problem. This approach is also suitable to compute higher order corrections of the 
dispersion law. Thus, we further collect 0(k) terms: 

r« = 2i\ ((1 - z) 1 - 1 /" - (1 - zj) . 

This function is integrable only for n > 1/2. Next, we compute p^: 



p „(1) rl 

(4-15) «* r = - T V-i I r^dz. 



(2) 

From this equation, we deduce that Txz is given by 

t$ = -i [ - r« + R e (U^ - (A<°> + U )(^y - AW(0(°))')) dz, 

(2) = -('(' - y)) 1/B -VS) + t^dydz. 

Jo Jo n 

We deduce that A^ 2 ) = (2) | 2= i is given by 

A ( 2) = R 2 ( 4 A^ 2 (8n 2 + 25n + 12) 

(2n + l) 2 (3n + 2) (4n + 3) (5n + 4) 

2 A^" 1 (5n 2 + 14n + 6) 



-ReCOt(6) 

(4.16) 4 



(2n+ l) 2 (3n + 2) (4n + 3) 
4 A" (n + 2) 



3 (2n- 1) (3ra+l) 
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To the next order one finds: 



r£> = -2i(A(l - z)) 1 - 1 ^)', 



F 2 



i? e (f/ o + A (o) )0 (o) . 



This yields the following expressions for (ri^, 0( 3 )) 
- 1 'iLp< 2 > 



-(3) 



F 2 



r£> + R e (U^ - (A® + C/ o )(0 (2) )' - A (1) (0 (1) )' - A (2) (0 {O) )O dz, 



{3) = - f f -(A(l - y)) 1 /"" 1 ^) + <jF>dydz. 
Jo Jo n 

As a result, we deduce that A^ = 0( 3 )| 2=1 reads 



K^ = Rl 



i % A^" 3 (408 n 5 + 3414 n 4 + 8935 n 3 + 9890 n 2 + 4905 n + 900) 
(2n + l) 3 (3n + 2) 2 (4n + 3) (5n + 4) (6n + 5) (7n + 6) 



-i^cot(0) 



+R e cot 2 (6) 



3 i A n- 2 (244 n 5 + 2248 n 4 + 5675 n 3 + 5940 n 2 + 2776 n + 480) 
(2n + l) 3 (3n + 2) 2 (4n + 3) (5n + 4) (6n + 5) 



2 i A"" 



(2n+ If (3n + 2) 



- cot(0) 



iAn (2n + 7) 
3 (2 n- 1) (4n+ 1) 



— iL 



iAn" 1 (988n 6 + 5666n 5 + 11622n 4 + 10825n 3 + 5008n 2 + 1120ra + 96) 
(2n - l)(2n + l) 2 (3n + l)(3n + 2)(4n + l)(5n + 2)(5n + 3) 



As a by product of this analysis, one recovers the classical values found in the 



expansion of the dispersion relation for a Newtonian fluid 



:K],p](n = 1, A = 3 if 



the reference velocity is the mean value of the Nusselt solution): 
A<°> = -3, 

A( 2 ) = -y J R e cot(0) + y J R e 2 + 3, 



A (3) 



. / 75872 



17363 



■ p3 _ 

V 25025 e 5775 



1413 



9 



R 2 e cot(6) + -R e cot 2 (6) + ^-R e cot(6) - Jcot(0) 
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4.4 Comparison with other shallow water models 

The shallow water models proposed in this paper are consistent up to order 2. In 
INM], a first order shallow water model was proposed but also shown inconsistent 
since it did not predicted the correct stability threshold (derivedfrom a consistent 
Benney type equation). Since then, other models was proposed in fRQ]. In particular, 
in the shear thickening case n > 1 where the regularization process does not introduce 
singularities, the authors proposed the shallow water model: 

d t h + d x q = 0, 

R e d t q = R e (-F(n)ld x q + G(n)^d x h) + I(n) ( h(l - cotm(9)d x h) - 



\n-i I t / ^?/ fl n 2 „ ,__sd x qd x h T ,^q 



+(^2) - K x {n)^±- - L 1 (n)^-d xx h + M 1 (n)d xx q ) , 



with 



n . lln + 6 . 6n + 3 T . . 3n + 2 

^0 = 7; — rT> = 1 — r^T> G ( n ) = 1 — r^7> H n ) 



2n + V w 4n + 3' v 1 4n + 3' v 7 2(2ra + l)' 
_ (2ra + l)(3n + 2)(12n 3 + 36n 2 + n - 1) (2n + l)(2n + 7)(3ra + 7) 

^ ~ 6n(2n- l)(3n+ l)(4n+ 1) ' ^ ~ 6(4n + l)(2n - 1) 

The coefficients i^i(n) and J\{n) are unimportant here since we are going to study the 
stability of steady states under small wavenumber perturbations. For that purpose, 
we linearize the previous system about (h,q) = (l,0o)- The system now reads 

d t h + d x q = 0, 

R e d t q = R e (G(n)(f)ld x h - F{n)(j)d x q) 

+I(n) {{2n + l)(h-q)- cotan(6)d x h) 
+M 1 (n)d xx q - L^fod^h. 

We search for plane wave solutions in the form (h,q) = e lk (H k ) t + x ) (}i(k),q(k)^ and 
expand both eigenvalues A and (h, q) with respect to k — > as 

2 

(A, h,q) = J2 k j {A U) , h {j \ q {j) ) + 0(k 3 ), 
3=0 
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with the normalization = 1 and = 0,j = 1,2. Inserting this expansion into 
the linearized system yields 

Ah + q = 0, 

iR e kAq = ikR e (^(n)<plh — F(n)(j)oq^j 

+I(n) ((2n + l)(h-q)- ikcot&n(9)hj 
+k 2 (L^foh - M 1 (n)qj . 

By collecting O(k ) terms into this system, one obtains 

A(°) + g ~(°) = 0, 

So that, one finds = —1. Next, by collecting 0(k) terms into this system, one 
finds 

AW+gM = 0, iRekWqW = iR e (G{n)<j>l - F{n)<j> ) -J(ra) ((2n + l)g (1) + icotan(0)) . 
From this system, we deduce 

A« = i C ° tan(g) - iR e7 % -. 

2n+l (2n + l)(3n + 2) 

At this step, one can conclude that the shallow water system found in ll^Q] is con- 
sistent up to order one (at least at the linearized level). At the second order, Ruyer 
Quil et al. pointed out that the terms quadratic in d x h, d x q found as n = 1 was not 
consistent with the order two shallow water model. Next, we show that even at the 
linearized level this system is not consistent except in the case n — 1. Indeed, by 
collecting 0(k 2 ) terms in the linearized system yields 

A(2) + gm = o, 

iR e (A^qW + A(o)g(i)) = -iR e F(n)(f) qW _ ( 2n + l)/(n)g (2) + L x {n)4> Q - M x (n). 
From this equation, we find 



A ( 2) = R 2 4(5n 2 + 14n + 6) 

+ 2) 2 (4^ 
2 (5n 2 + 14n + 6) 



(2n+ l) 2 (3n + 2) 2 (An + 3) 



R e cotan(6)- 

(2n + l) 2 (3n + 2) (4n + 3) 

4 (n + 2) 



3 (2n- 1) (3n + l)' 
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One finds indeed that, in comparison to the A^ 2 ^ of the Orr-Sommerfeld computations, 

(2 n + 1) (3 n + 2) (4 n + 3) (5 n + 4) 

As a consequence, the model derived in fSoj for shear thickening fluids (n > 1) not 
consistent, even at the linear level with Navier-Stokes equations.: as quoted in KQ], 
this is due to the fact that some corrective terms due to inertia are neglected. 

5 Conclusion and Perspectives 

In this paper, we have derived consistent shallow water equations for thin films of 
power law fluids down an incline. For that purpose, we have introduced a new, 
weaker, formulation of the Cauchy momentum equations. In particular, we have 
kept the fluid strain rate as an unknown in order to deal with the divergence of the 
apparent viscosity near the free surface. In this new framework, we have computed 
an expansion of the velocity field, in particular the bottom strain rate and streamwise 
discharge rate up to order 2 with respect to e, the film ratio. We have then generated 
a family of consistent models of shallow water type. Up to first order, we proposed a 
new model that is consistent, conservative and invariant by Galilean transformation. 
Up to second order, we proposed the first fully consistent shallow water model that 
is an extension of the one proposed in |RQ] . 

One classical test to verify the consistency of shallow water models is to compare 
the spectrum associated to the linearized equations about a constant state and the one 
that arises from a direct analysis of the Cauchy momentum equations. As a byproduct 
of our new formulation of Cauchy momentum equations, we propose the first rigorous 
spectral analysis of power law fluids down an incline. We have derived a suitable 
linear problem that is an generalization of Orr Sommerfeld equations for newtonian 
fluids. We have shown that the nonlinear expansion of the velocity field in the shallow 
water scaling provides a natural expansion of eigenvalues and eigenfunctions of this 
generalized Orr Sommerfeld equations. 

Our analysis is restricted to n > 1/2 and strongly suggests that the free surface 
Cauchy momentum equations are ill posed if n < 1/2 even at the linearized level. As 
a consequence, the use of the power law in some applications like glaciology n — 1/3 
or for various polymers is ques ionable: in such situations, one has to regularize the 
power law at low strain rate. In KQ], the authors used such a regularization procedure 
in order to compute shallow water models for shear thinning fluids n < 1 whereas 
it is useless if n > 1. Here, we show that both cases are treated similarly provided 
that n > 1/2: as a result, one has to recover our models and expansions when the 
parameter of regularization goes to zero if n > 1/2. In contrast, if n < 1/2, one needs 



REFERENCES 



33 



a regularization procedure to derive shallow water models and these models should 
strongly depend on the regularization procedure: this question shall be treated soon 
in a future work. 

There are various exemples of regularization of the power law founc in the liter- 
ature like the Carreau law or mixed power law/newtonian law used in KQ]. In the 
later case, a strain threshold is introduced: below this threshold, the fluid is newto- 
nian and above, it is a power law. The authors introduced this law in order to do 
computations explicitly but the threshold is fixed which may leads to some inconsis- 
tency if one let this threshold goes to zero. The Carreau law seems more natural but 
computations are not explicit. Then, a generalization of this law that would lead to 
explicit computations is of interest in order to compare with our results: this problem 
and the influence on the derivation of shallow water models, stability properties of 
constant states are very interesting for applications purposes (sea of ice, flows of poly- 
mers with n < 1/2). Another interesting extension to this work would be to study 
the influence of the boundary conditions at the bottom, in particular if the no-slip 
condition is replaced by the Navier-slip condition: this is of interest if one considers 



Uiipn: 
1UMBR 

flows of power law fluids down porous media IfOTvlBR] 
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